Cohesive energy, stability and structural transitions in polyelectrolyte bundles 
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A lattice of uniformly charged, infinitesimally thin, rods decorated with an ordered array of 
counterions exhibits anomalous behavior as the spacing between the rods is varied. In particular, the 
counterion lattice undergoes a sequence of structural shearing, or "tilting," phase transformations 
as the spacing between the rods decreases. The potential implications of this behavior with respect 
to the packaging of biologically relevant polyelectrolytic molecules are commented upon. 

PACS numbers: 61.20.Qg, 61.25.Hq, 64.70.Kb, 77.80. Bh 



I. INTRODUCTION 



Polyelectrolyte chains naturally repel each other since 
they carry an overall electric charge. That these chains 
will nevertheless form condensed phases in the presence 
of oppositely charged counterions has been known for 
some time. It has also become clear that this condensa- 
tion results from some form of organization of the coun- 
terions, either dynamical [3, I2| essentially static, in 
the form of a counterion lattice In both cases, the 
attraction between oppositely-charged regions overcomes 
Coulomb repulsion between the bare chains and leads to 
a net attraction. Correlations in the counterion system 
are crucial, since it is well-known that mean field the- 
ory, represented by solution of the Poisson-Boltzmann 
system of equations, cannot yield an attraction between 
like-charged rods dH H ■ 

The physics of polyelectrolytes is relevant to biological 
systems. For example, double helix DNA will, under cer- 
tain conditions, organize into a condensed state in which 
it self-assembles into bundles of densely packed parallel 
rods U, [31 • An important venue for this condensation is 
within the heads (or capsids) of various viruses 1^. This 
organization occurs in spite of the fact that "naked" DNA 
carries a strong negative charge — one excess electron per 
phsophate group on the backbone . DNA has also been 
observed to form condensed liquid-crystal-like phases in 
the presence of polyvalent counterions XQ_^ .lli] . In ad- 
dition, counterion-mediated formation of actin bundles 
has also been observed Interestingly, in this last 

case, there is evidence of a kind of counterion lattice, in 
the form of a charge-density-wave-likc modulation along 
the axis of the condensed actin filaments. In such a lat- 
tice, the "counterions" in the lattice actually consist of 
clusters of individual counterions, and can therefore not 
be represented as point particles. The question of an 
ionic lattice specifically in the case of condensed DNA 
remains an open one, since there is, as yet, no experi- 
mental indication of such an organization of counterions. 
Furthermore, theoretical estimates indicate that a 3D 
structure requires counterion valences, Z, greater than 



about 6 1^]. A lattice of point-like counterions in con- 
densed biological rodlike molecules thus remains a con- 
jecture, rather than an established fact. Nevertheless, it 
represents an sufficiently interesting and potentially im- 
portant construct within biological and polyelectrolyte 
physics that its properties merit investigation. 

In this paper we analyze more carefully the cohesive en- 
ergy of the background of negatively charged rods along 
with a neutralizing, crystalline counterion system. We 
assume the counterions form a modified face-centered 
cubic (fee) crystal. Our principal result is that, as the 
rod density increases, the counterion system maintains 
stability by undergoing structural transitions to "tilted" 
lattices. More precisely, we have located two symmetry- 
breaking structural transitions: the first is of the three- 
state Potts type and is therefore weakly first order. 
The second is a continuous transition in the Ising uni- 
versality class. A fiuctuation analysis indicates that the 
transitions are mean field in nature. This reflects the 
dominating influence of long-range Coulomb interactions 
between the unscreened charges that are bound to the 
rods. At the end we speculate on how the singulari- 
ties associated with such transitions could be relevant 
to the physics of packaging and other bundling phenom- 
ena in biology. It is worth noting that structural tran- 
sitions have been recently observed in DNA-dendrimer 
complexes 0]. These transitions do not appear to cor- 
respond to the speciflc ones we discuss. However, the 
notion of structural transitions in complexes of rods and 
localized charges clearly has an experimental as well as a 
theoretical basis. 

We also touch on, but do not discuss in detail, the 
counterion "melting" transition expected when the rods 
are far apart. This leads to a "counterion liquid." The 
attraction between rods in such a state has been exten- 
sively explored by Ha and Liu, 0, ITsl IT^ . The melting 
transition is continuous and can be shown to be in the 
universality class of the three-dimensional XY model; 
see, e.g., [13 111, IS ■ A brief review of the arguments 
leading to this conclusion appears in an appendix. 

An outline of this paper is as follows. The following 
section contains the characterization of the lattice and 
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the Coulomb energy calculations. Section Hill addresses 
"melting" of the counterion lattice, while Section llVl con- 
tains concluding remarks. A series of appendices address 
some technical issues. Appendix A contains a calculation 
of the Coulomb energy of a lattice of infinite, uniformly 
charged rods. Appendix B provides the derivation of the 
term in the energy quadratic in displacements of coun- 
terions from their lattice sites (i.e., harmonic "phonon" 
dispersion relation). Appendix C provides some heuris- 
tics on the "melting" of the counterion lattice. 



II. LATTICE OF RODS AND COUNTERION 
ENERGY 



The "lattice" consisting of a set of negatively charged 
rods and attached polyvalent counterions of charge ~\-Ze 
will be treated as a simple ionic lattice. At the first stage 
of approximation, the fixed charge of the rods will be re- 
placed by a uniform negatively charged back grou nd (as 
in the jellium model of interacting electrons |21|). The 
counterions will be assumed to form a three-dimensional 
lattice that generalizes the close-packing arrangement 
in a face-centered cubic (fee) crystal. Given that the 
rods are infinitesimally thin, there is no prospect of a 
two-dimensional arrangement of charges on the surface 
of any of them As has been noted in the litera- 
ture [23 , such a lattice can be constructed starting from 
a two-dimensional hexagonal close-packed structure, via 
the introduction of three sublattices. Our generalization 
makes use of this construction. No significant error is 
introduced by our neglect of the explicit contribution of 
the rods themselves to electrostatic interactions. Dimen- 
sional considerations and detailed calculations (see Ap- 
pendix^ lead to the conclusion that the Coulomb en- 
ergy due to interactions between charged rods consists 
of two contributions, the first independent of the separa- 
tion between rods and the second going as the logarithm 
of the separation between rods, assuming overall charge 
neutrality. Given that the counterions are forced to sit 
on the rods, the principal outcome of the rod-counterion 
interaction arises from the overall charge neutrality en- 
forced by the charge on the rods. 



A. Structure of the rod-counterion lattice 

Imagine a large bundle of hexagonally close-packed 
pencils or rods. Viewed end-on, the rods lie on three 
triangular sublattices, as depicted in Fig. ^ Now pass 
planes at regular spacing ly perpendicular (initially) to 
the rods. Where the first plane intersects the rods on the 
first sublattice, place counterion charges. Likewise treat 
the second and third planes at their intersections with 
rods on the second and third sublattices, respectively. 
Repeat this pattern with subsequent planes, thereby 
building an infinite three-dimensional modified fee lat- 
tice. Figure 121 displays the lattice from a viewpoint per- 




FIG. 1: Schematic of the counterion lattice, seen end-on with 
respect to the rods that support it. The rods are hexago- 
nally close-packed, and the three sublattices occupied by the 
counterions are indicated. 



pendicular to the axis of the rods. The three counterion 
sublattices are evident if one scans the figure from left to 
right. If Z„ is adjusted appropriately relative to the lattice 
constant of one triangular sublattice, lh^ the counterion 
charges themselves sit on a true fee lattice. 




FIG. 2: A view of the lattice perpendicular to the axis of the 
rods. The spheres represent the counterions. Scanning from 
left to right, one observes the three sublattices. 

The three primitive vectors of the Bravais lattice of 
counterions so constructed are 



a\ = 

3.2 = 

a-s = 



^IhX + + Iv (^9x + gyVs^ z (2.1) 

-^kx + ^IhV + lv (^-(Jx + gyVsj z (2.2) 

(2.3) 



The dimensionless vector g ~ xgx + ygy encodes the 
possibility of tilting the planes of counterions. The cor- 
responding primitive vectors of the reciprocal lattice, 6i, 
are constructed in the standard way |23|, so that the re- 
lationship between the two sets of primitive vectors is 
di ■ bj = l-nbiA ■ To recover the true fee lattice, the aspect 
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ratio, r = Ih/lv, is taken to be 1/8/2 and g is set equal 
to zero. The aspect ratio will play the role of a control 
parameter in what follows. 

Assume overall charge neutrality of the rod-counterion 
system, so that there is exact cancellation between the 
mean charge per unit volume of the counterion lattice and 
the uniformly charged negative background provided by 
the rods. We make use of the Ewald method for the eval- 
uation of the Coulomb sum [Ulillil. The "reduced" 
Coulomb energy is defined as 



= -£'Coul 

= E 



S (2.4) 



where n is a triplet of integers (711,712,713), and the sub- 
traction S represents the compensating interaction with 
the smeared out negatively charged background. 

In our calculation, we assume that the distance be- 
tween "planes" of ions is fixed at Z„. The effect of a 
compression of the lattice of polyelectrolytic rods is to 
decrease the spacing, Ih- 



B. The Coulomb sum 

Using Ewald summation techniques |^ |22, H^, one 
generates an expression for the Coulomb energy that can 
be expressed in terms of the sum of four terms. Those 
terms are 



-1/3 



(2.5) 



E 



\niai + n2a2 + 7130.3] 



X (^1 - erf (^iTiifli -f- 712^2 + "-303 {llQ ^^'^) ) 

(2.6) 



—a 

V 



2/3 



(2.7) 



4_- „-\mibi+rn2b2+m3b3\^a-^ {lll^y^yi 
— ^ (2.8) 



\mibi + 771262 + 'Tis&sp 



The quantity a in the expressions above is an adjustable 
parameter, which is ideally set equal to a value that max- 
imizes convergence of the Ewald sums in (|2.6|l and H2.8|l . 
A close-to-optimal choice is a = An. 

The mean-field phase diagram of the lattice can be 
determined by examining the dependence of the Coulomb 
energy on the aspect ratio, r, and the tilt vector, g. What 
one finds is that for r greater than a threshold value w 



1.1 the Coulomb energy is minimized when g = 0. At this 
threshold value, three minima lying symmetrically in the 
5-plane represent equally low energies. This situation is 
illustrated in Figure |31 When r < 7^, g = no longer 



0.04 
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FIG. 3: Contour plot of the Coulomb energy of the generalized 
fee counterion lattice in the — Qy plane at r — ~ 1.1. 
The four minima, including the one at g = 0, are all of equal 
depth. 

represents a global minimum of the energy. At r = « 
1.097, the local minimum at g = disappears. In this 
sense, one can think of "spinodal" point. Note 

the small difference between Ts and Va ■ The transition at 
Ta is weakly first order. 

As r is further reduced, corresponding to even closer 
packing of the rods, a new structural transition is en- 
countered, at which the three minima each split into two 
new ones. The aspect ratio, r^, at which this transi- 
tion takes place is approximately equal to 0.801. Con- 
tour plots illustrating the onset of this transition and 
the evolution of the new minimum energy configurations 
arc shown in Figure 01 That this transition is also a 
symmetry-breaking one is evident from Fig. 0] Here, the 
transition is continuous, in that the appearance of the six 
minima at lower aspect ratio is simultaneous with the dis- 
appearance of the three minima associated with values of 
r greater than rf,. At mean field level this transition is of 
the standard second-order type. Given the nature of the 
symmetry-breaking, it is altogether reasonable to classify 
it as an Ising-like, or 0(1), phase transition. In the case 
at hand, the aspect ratio plays the role of temperature. 
This means that singularities that one expects to find in 
the entropy and specific heat of a thermal system will 
show up here in the form of nonanalyticities in the de- 
pendence of the energy as a function of spacing between 
rods, with direct consequences on packing forces. 

In light of the continuous nature of the phase transition 
at r = Tf, and the weakness of the first-order phase tran- 
sition at r — ra, the question of the effects of ffuctuations 
is clearly relevant. We address this question with the use 
of coarse-grained effective Hamilitonians. Assume that 
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FIG. 4: The emergence of six minima in the Coulomb energy, 
and the migration of those minima as the aspect ratio is de- 
creased below the threshold value ri, = 0.801. Reading from 
left to right and top to bottom, the values of r are equal to 
0.801, 0.693, 0.577 and 0.433, respectively. The arrow in the 
bottom-right-hand figure points to one of the six true minima 
of the energy. 



the counterions are only allowed to move along the rods. 
Then, at the onset of the first-order transition, we can 
express the fluctuations in the locations of the counte- 
rions in terms of a scalar displacement field, u{r) that 
quantifies counterions displacments paralles to the rods. 
The Coulomb energy can then be written as an expansion 
in terms of u. This energy is most usefully expressed in 
terms of the Fourier transform of the displacement field. 
The most relevant terms in the expansion yield the ex- 
pression 



H = 



k\^{r-rs) + Ck^ + AQ^ + B 



E 

k,Q 

xu{k, Q)u{—k, — Q) 



fc2-HQ2 



E ^3(<fl, 92, g3)"(gi)u(g2)w((f3)%+g2+«3 

.■■■,g3 

E "'4(91, 92, 93, q4)u{qi)u{q2)u{q3)u{q3) 

(2.9) 



gi,...,g4 
xSff, 



91+92+93+94 

Here, we have split the three-dimensional wave vector q 
into a two-dimensional vector, k in the x-y plane and a 
component, Q, in the z-direction. The term with coef- 
ficient B is the contribution to the quadratic energy re- 
fiecting the long-range nature of the unscreened Coulomb 
interactions between counterions. We have inserted the 
"stabilizing" quartic term going as fc*, but have ignored 
the inessential term proportional to Q^, or the cross term. 



proportional to k^Q^. It is important to keep in mind 
that the actual stabilization of the system results from 
the shearing transition. 

The higher order terms in the the effective Hamiltonian 
(|2.9() have the following forms: 



and 



"'3(9*1, 92, 93) = -Ci (kl 



Wiiqi, • • • ,94) 

= W{qi,---,q4) 



3kyk^j 



(2.10) 



{kl ■ k2){k3 • ^4) -t- permutations 



+B2Q1Q2Q3Q4 



+B3 



(kl ■ ^2)^3(94 + permutations . (2.11) 



Here, 



W{qi,---,q4) 

= viqi + q2) + v{qi + 9*3) + v{q2 + 93) 

-v{qi) - v{q2) - viqs) - v{-qi - 92 - 93)(2-12) 



where 



v{q) 



91 

I9T 



(2.13) 



There are three distinct permuations of the indices in the 
first line of the right-hand-side of Eq. (|2.11|l and six dis- 
tinct permutations in the last line of the right-hand-side 
of that equation. The term Vl^(9i, . . . , qi) appears to be 
the most relevant contribution to the fourth-order cou- 
pling in the Ginzburg-Landau- Wilson model appropriate 
to this system. However, as it turns out, the most im- 
portant term in Eq. (|2.11() is the first one in the square 
brackets, going as fc^. 

That the contributions to the energy associated with 
fluctuations about the lattice have the forms shown 
above can be established through explicit evaluation of 
Coulomb-type lattice sums. Appendix^outlines the cal- 
culation in the case of of the quadratic terms and presents 
results for that and the fourth order term that are ob- 
tained through explicit evaluation of those sums. 

To perform an analysis of the effects of fluctuations 
on this weak first-order transition, one can consider a 
"Ginzburg criterion" ^ applied to the one-loop contri- 
bution to the "entropy" of the system [2^. Recall that 
the leading contribution to the mean field entropy goes as 
r — Ts for r — s- rs + . Given the form of the of the quadratic 
term in Eq. 1)2. 9|l . and recalling that s ^ dF/dr, with F 
the free energy, we have the following expression for the 
one-loop entropy. 



k^d^ 



q 



[r- 



rs)k^ + Ck'^ + AQ"^ 
Pd^kdQ 



{r --r,)k^ + Ck^ + B 



(2.14) 
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The expression on the right-hand-side of Eq. (|2.14() con- 
tains the important terms in the denominator. A variety 
of methods exist for the evaluation of the integrals in 
this expression. The essence of the results follows from a 
rescaling Q = k^x. Then, the integral to perform is 



k'^tPkdx 



(r - rs) + Cfc2 -f Bx"^ 



(2.15) 



E3 
O 

u 




1.59 



We now note that this has the same qualitative depen- 
dence on the "reduced temperature," r — Ts, as the cor- 
responding one- loop integral of a five-dimensional 0{n) 
model with short range interactions. A further rescaling 
of the integration variables produces a leading singularity 
proportional to (r — r^)'^/^ compared to the mean- field re- 
sult oc r — . This application of the Ginzburg criterion 
shows that the transition is unrenormalized. A calcula- 
tion of the renormalization of the fourth-order interaction 
leads to the same conclusion, namely that fluctuations 
lead to a well-behaved change in the amplitude. This 
leads to the conclusion that, because of the long range 
dipole-dipole interactions between the charges embedded 
on the rods, the first-order shearing phase transition of 
the charge lattice at r = is essentially mean field in 
nature. Similar arguments, allowing for modifications in 
the form of Eq. H2.9|) by, for example, breaking rotational 
symmetry in the x — y plane, reveal that fluctuations 
don't modify the mean-field continuous shearing transi- 
tion. Technically, the transition at r = r;, is Ising-like 
with effective dimensionality d > 4. 

The full results of the Coulomb energy calculations are 
shown in Fig. [S] The analysis above reveals that fluctu- 
ations don't qualitatively change the picture. 



Coul 




FIG. 5: Solid line: the reduced Coulomb energy as a function 
of aspect ratio in the vicinity of the first-order and continuous 
transitions. Dashed line: the energy if the counterion lattice 
is constrained not to tilt. The locations of the transitions are 
indicated on the figure. 



Figure 1^ displays the (very small) discontinuity in the 
derivative dEcoui/df at the first order transition at r = 



FIG. 6: A plot of dEcoiii/dr in the immediate vicinity of 
r — ra, illustrating the discontinuity there, associated with 
the first order phase transition to the sheared lattice. 



III. "MELTING" OF THE COUNTERION 
LATTICE 

The three-dimensional lattice of statically correlated 
counterions being discussed can be reasonably expected 
to exist at low enough temperatures and for sufficiently 
close packing of the lattice of rods. At high temperatures 
and when the spacing between rods is relatively large, 
the counterion lattice will not be stable against thermal 
fiuctuations. A variety of arguments lead to criteria for 
the existence of this lattice. The simplest compares the 
(negative) interaction energy of the lattice to ksT. It 
is reasonable to expect that the lattice will resist the 
disordering effects of thermal fiuctuations when 



Coul 



(3.1) 



Assuming that the general condition, Z'^e'^ /{ksTly) ^ 1 
for the existence of a "Wigner crystal" in the couterion- 
polyelectrolyte system are met > expect the crystal 
discussed here to be stable as long as l-Ecouil is of order 
unity or greater. According to Fig. [3 this should be 
the case as for a range of aspect ratios that encompass 
the two structural transitions on which this paper has 
focussed. 

An alternative estimation of the threshold at which 
melting of the counterion lattice takes place is based 
on the Lindemann criterion [JJ, according to which a 
lattice melts when thermally-induced displacements are 
some fraction of the lattice spacing. This leads to ap- 
proximate predictions for the melting transition that we 
expect to be consistent with the energy- us.-entropy ar- 
guments leading to the criterion (|3.1|) for the stability of 
the counterion lattice against melting. 

When melting takes place, it does so in the same way 
that charge density waves disappear, that is to say, con- 
tinuously, with thermodynamic signatures that identify 
its universality class as that of the 3DXY model. The- 
oretical arguments and experimental observations that 
justify this conclusion in the case of charge density wave 
systems can be found in the literature [l^ [T3. [l9l l20j . 
Appendix El contains a brief, heuristic argument for the 
nature of the melting transition in the counterion system. 
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IV. CONCLUSIONS 

As noted above, the singularities in the dependence of 
the Coulomb energy on aspect ratio have implications 
with respect to the energetics of bundling or packag- 
ing polyelectrolytic molecules. These implications fol- 
low from the effective force exerted by the decorated 
rods on each other, as determined by the "pressure" , 
p = —dEcoui/dr. The pressure (analogous to the en- 
tropy in a thermal system as discussed) has a disconti- 
nuity at r = Ta, and a discontinuous derivative at r = r^. 
The singularities are not striking (recall Fig. EJ. The 
kinetics of bundling or packaging such molecules with 
crystalized counterions are presumably affected by the 
Coulomb energy, and should, in principle, reflect these 
singuarities. What is more important, however, is the 
possibility that a sequence of structural transitions in the 
counterion lattice assists compression to high densities. 
In Fig. 13 we have shown the Coulomb energy computed 
directly from the Ewald sums in the vicinity of the first- 
order transition at r = and the continuous transition 
at r = rb- The dashed curve shows the energy if the 
lattice were constrained not to tilt. It is interesting to 
contemplate whether in some instances Nature relieves 
the strong Coulomb repulsion via structural transitions. 

It is important to note that one cannot argue for struc- 
tural transitions as the sine qua non of polyelectrolyte 
packaging. Specifically, in the packaging of DNA in tail- 
type bacteriophages, ATP is known to provide the fuel for 
a packaging motor ^ . Order of magnitude estimates 
suggest that ~ 50 — dQpN forces generated by this mo- 
tor |2_6] suffice to overcome Coulomb forces and compress 
a rod-like system to observed densities without benefit of 



structural transitions. Of course, in experimental situa- 
tions other repulsive energies besides Coulomb are also 
involved. 

The present analysis does not suggest dramatic conse- 
quences in actin bundling, or in the kinetics of packaging 
and/or infection in a bacteriophage life-cycle. However, 
a continuous shearing transition, if realized, should be 
accompanied by strong counterion charge fluctuations, 
which could be susceptible to dynamic scattering exper- 
iments. The tilting, estimated to be a few percent at 
the first-order transition, could potentially be detected 
by standard diffraction techniques. 

While, strictly speaking, we have shown that a shear- 
ing transition ought to occur in a sufficiently closely- 
packed, constrained polyelectrolytic system, we have not 
ruled out the possibility that some other transition inter- 
cedes, preempting this particular rearrangement as the 
system of rods is compressed. One would have to con- 
sider possibilities for ever larger unit cells in the rod lat- 
tice. More importantly, we can only suggest the possibil- 
ity that tilting or other structural transition(s) occur as 
biologically relevant polyelectrolytic molecules condense. 
We hope this work will stimulate further experiments on, 
for example, the kinetics of DNA packaging, since such a 
transition could occur "on the fly" as the "spooling" pro- 
gresses. i2^^] In other potential experiments, DNA or 
other bundles of varying density could be prepared and 
probed statically. 

We thank Prof. A. B. Harris, Dr. Ron Fisch and 
Prof. Roger Hendrix for helpful comments. We are also 
grateful to Prof. Robijn Bruinsma for bringing Ref. [T^ 
and the work on which it reports to our attention. 



[1] F. Oosawa, Polyelectrolytes (M. Dekker, New York,, 
1971). 

[2] B. Y. Ha and A. J. Liu, Phys. Rev. Lett. 79, 1289 (1997). 

[3] B. I. Shklovskii, Phys. Rev. Lett. 82, 3268 (1999). 

[4] N. Gronbech- Jensen, R. J. Mashl, R. F. Bruinsma, and 

W. M. Gelbart, Phys. Rev. Lett. 78, 2477 (1997). 
[5] F. Oosawa, Biopolymers 6, 134 (1968). 
[6] T. Ohnishi, N. Imai, and F. Oosawa, J. Phys. Soc. Jpn 

15, 896 (1960). 
[7] V. Bloomfield, Biopolymers 31, 1471 (1991). 
[8] R. W. Hendrix, Cell 94, 147 (1998). 
[9] J. D. Watson, Molecular biology of the gene (Ben- 

jamin/Cummings, Menlo Park, 1988), 4th ed. 
[10] D. Durand, J. Doucet, and F. Livolant, Journal de 

Physique II (Atomic, Molecular and Cluster Physics, 

Chemical Physics, Mechanics and Hydrodynamics) 2, 

1769 (1992). 

[11] J. Pelta, D. Durand, J. Doucet, and F. Livolant, Bio- 
physical Journal 71, 48 (1996). 

[12] T. E. Angelini, H. Liang, W. Wriggers, and G. C. Wong, 
Proc Natl Acad Sci U S A (2003). 

[13] D. J. Amit and A. Shcherbakov, J. Phys. C (Solid State) 
7, L96 (1974). 



[14] H. M. Evans, A. Ahmad, K. Ewert, T. Pfohl, A. Martin- 
Herranz, R. F. Bruinsma, and C. R. Safinya, Struc- 
tural polymorphism of DNA-dendrimer complexes, un- 
pubhshed. 

[15] B. Y. Ha and A. J. Liu, Physical Review E (Statistical 
Physics, Plasmas, Fluids, and Related Interdisciphnary 
Topics) 58, 6281 (1998). 

[16] B. Y. Ha and A. J. Liu, Physical Review E (Statistical 
Physics, Plasmas, Fluids, and Related Interdisciphnary 
Topics) 60, 803 (1999). 

[17] G. Griiner, Density waves vn solids, Frontiers in physics 
; v. 89. (Addison- Wesley, Reading, 1994). 

[18] S. Girauh, A. H. Moudden, and J. P. Pouget, Phys. Rev. 
B 39, 4430 (1989). 

[19] H. Fukuyama and H. Takayama, in Electronic properties 
of inorganic quasi-one-dimensional compounds. Part I: 
Theoretical, edited by P. Monceau (D. Reidel ; Sold and 
distributed in the U.S.A. and Canada by Kluwer Aca- 
demic Publishers, Dordrecht Holland ; Boston Hingham, 
MA, U.S.A., 1985), pp. 41-104. 

[20] J. W. Brill, M. Chung, Y. K. Kuo, X. Zhan, E. Figueroa, 
and G. Mozurkewich, Physical Review Letters 74, 1182 
(1995). 



7 



[21] J. M. Ziman, Principles of the theory of solids (Cam- 
bridge U. P., New York, 1979), 1st ed. 

[22] N. W. Ashcroft and N. D. Mermin, Solid state physics 
(Holt, New York, 1976). 

[23] P. P. Ewald, Ann. Physik 64, 253 (1921). 

[24] V. L. Ginzburg, Sov. Phys. Solid State 2, 1824 (1961). 

[25] Y. M. Ivanchenko and A. A. Lisyansky, Physics of criti- 
cal fluctuations, Graduate texts in contemporary physics 
(Springer- Verlag, New York, 1995). 

[26] D. E. Smith, S. J. Tans, S. B. Smith, S. Grimes, D. L. 
Anderson, and C. Bustamante, Nature 413, 748 (2001). 

[27] S. C. Riemer and V. Bloomfield, Biopolymers 17, 785 
(1978). 

[28] T. Odijk, Biophys J 75, 1223 (1998). 
[29] L. W. Black, in The Bacteriophages, edited by R. Calen- 
der (Plenum, New York, 1988), vol. 2, pp. 321-373. 



APPENDIX A: THE COULOMB SUM FOR 
UNIFORMLY CHARGED RODS 

As discussed in the text, the charges on the rods have 
been assumed to be smeared out into a three-dimensional 
smooth background, analogous to the jellium model used 
in discussions of electron-electron interactions in a metal 
PH . In this Appendix, we evaluate the corrections to 
this approximation by now assuming that the backbone 
charges are uniform along the rods. The method utilized 
to evaluate the alteration in the Coulomb energy that is 
induced by this refinement in the model is a version of the 
Ewald summation technique |^ |22, H^l . The quantity 
we will calculate is the potential energy of a charged test 
rod in the presence of an array of uniformly charged rods, 
which are assumed to be in a hexagonal close-packed ar- 
rangement. At the end, the test rod will be moved onto a 
rod of the lattice. Before doing this, we eliminate ("sub- 
tract") the interaction between the test rod charge and 
the rod on which it eventually sits. 

Because the rods are uniformly charged in the z- 
direction, the potential at the two-dimensional location 
f due to a rod at the origin is given by 



1 



(r) oc — 

TT 



iq-r -| 
TT 



dt 



(Al) 

We split the f-integration into one from to T and an- 
other from T to oo, where T is, initially, arbitrary. We 
suppose that the distance between nearest-neighbor rods 
in the close-packed lattice is Ih- Given this, we will take 
T (X The primitive vectors for this lattice can be 
taken to be 



0.2 



-y, 



2 2 

and the volume of the primitive cell for this lattice is 



(A2) 
(A3) 



vws 



V3 



(A4) 



The corresponding reciprocal lattice vectors are 



47r 



y, 



(A5) 
(A6) 



while the volume of the primitive cell of the reciprical 
lattice is given by 



VBZ 



87^ 



(A7) 



For the integration from t = T to t = oo, we make use of 
the two-dimensional version of the Poisson sum formula: 
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r=Ri 
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f{r)d'r, 



(A8) 



where {Ri] correspond to a hexagonal lattice, and where 
the Qk range over the reciprocal lattice. We then have 
for the contribution of this region of ^-integration to the 
potential energy of a point charge at the location r 
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d^qe 
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Qk 



g-TQi^iQ^-r 
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The only tricky part here is the term for which Q^. = 0, 
corresponding to the smeared-out portion of the distri- 
bution. This can be handled by a careful subtraction. 
Because of the fact that T cx l\ and the fact that the 



Qfc's go as 



91 are independent of the 



terms in 
spacing between rods, Ih- 

The next portion of the integration is from to T. 
Here, we sum directly in real space. We find 



E 

Ri 



1 t 



(AlO) 



A careful subtraction eliminates the divergent contribu- 
tion of term in the sum in which Ri ~ f. The divergence 
arises when the limit r — )■ is taken, corresponding to 
placing the test rod on one of the rods in the lattice. 
One subtracts the standard energy associated with the 
Coulomb interaction between the test rod charge and the 
rod at the origin. This subtraction leads to the only term 
with any 1^ dependence, going as InZ^. All other depen- 
dence disappears because of cancellations between the 
Z^-dependence of the Qk's, the i?/'s and T. 

The subtraction in (|A9|) associated with the Qk = 
term is easily carried out, as is the subtraction in the case 



8 



of the term associated with a test rod position r directly 
on a rod, in lAlOl) . Choosing specifically 



4^ ' 



(All) 



with a an arbitrary constant, one can perform the sums 
indicated above, with the appropriate subtractions. One 
finds the following result for the Coulomb sum: 



- 2.78608 + 2 In/,, 



(A12) 



energy associated with those distortions 



_ Q 

(B2) 

where u{q) is the spatial Fourier transform of the lattice 
displacement Uj: 



-iq-Rj 



(B3) 



Note that the above result must be independent of the 
arbitrary constant a. We are now in a position to include 
more precisely the effect of the concentration of negative 
charge on the polyelectrolyte rods. One must multiply 
the result (|A12p by 



Given that the interactions are Coulomb, the sum of 
interest is of the form 



(B4) 



Ze 



(A13) 



where a is the linear charge density on the rods, Z is 
the valence of the condensed counterions, and ly is the 
spacing between counterions on the rod. The above mul- 
tiplication will yield the energy of the rod lattice is per 
unit rod length. 

We are thus led to an expression for the contribution 
of the rod-rod interaction to the Coulomb energy; it has 
a non-trivial, but smooth dependence on the separation 
between rods. This energy adds as a "background" term 
to the energies calculated in the body of this paper. A 
very similar calculation is used to compute the interac- 
tion energy between the counterion lattice and the uni- 
formly charged rod lattice. Both contributions have no 
effect on the structural transitions discussed in the body. 



APPENDIX B: EWALD SUM FOR THE SECOND 
AND HIGHER ORDER TERMS IN 
COUNTERION DISPLACEMENTS 

The result for the energy cost of a distortion of the 
counterion lattice follows straightforwardly from the ex- 
pression for the energy of a collection of interacting 
charges. We start by writing the two-particle interaction 
energy in terms of its spatial Fourier transform. 



V{t) 



d^qv{q)e 



tq-r 



(Bl) 



Then the positions of the counterions, r*,-, are expanded 
in terms of displacements from the lattice sites, u{Rj). 
Assuming that those displacements are entirely in the z- 
direction, in line with our model in which the counterions 
are bound to the charged polyelectrolye rods, the displac- 
ment of the j'^ counterion from its equilibrium position 
on the counterion lattice will be equal to zu{Rj) = UjZ. 
Expanding to second order in the u's, we obtain for the 



where the Q's are displacement vectors on the reciprocal 
lattice. The primitive displacement vectors on this lattice 
are 6i, 62 and 63. The primitive displacement vectors on 
the original lattice are oi, a2 and 03. The relationship 
between the 5's and the a's is 



bl = 2tt 



a2 X 03 
ai ■ {0.2 X as) 



(B5) 



and similarly for 62 and 63. The primitive volume in the 
reciprocal lattice is the volume of the first Brillouin zone, 
given by 



vbz 



fei • (62 X 63 



(B6) 



The primitive volume in the real lattice is the volume of 
the Wigner-Seitz cell, given by 



vws = \ai ■ {0,2 X 03)1 



(B7) 



Given the relationship between the a's and the 6's, the 
the following holds 



vbzvws = (27r)"' 



(B8) 



The Poisson sum formula in three dimensions takes the 
following form: 
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iRi 



■fiQ) (B9) 



The sum on the right hand side is over all lattice points 
on the real space lattice. One final relationship between 
the a's and the b's is 



(BIO) 



In (|B9p each Q is of the form mi6i+m2fe2+"^3^3, where 
the mi's take integral values from —00 to 00. Similarly 
the i?'s are of the form riidi + 77.202 + n^a^, where the 
n^'s range over all integers, as well. 
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Now, one can clearly write the expression in (|B4|) in 
the form f{q) — f{0). Focus on f{q) and write 



(q. + QzY 



q + Q 



{qz+Q.fe-\^+^\"'dt (Bll) 



As in Appendix A, the integral over t splits into an inte- 
gral from to r and from T to oo, where, in this case, 
we choose 



T = /3 



An 



(B12) 



where (3 is arbitrary. For the integral from T to cx), one 



finds 



Q \q + Q 



2/3 
BZ 



E 



(B13) 



The Poisson sum formula is applied to the integration 
over t from to T. This leads to the following expression 



— T. 

VBZ ^ 
R 



2/3_^, 2/3 



^■R-\q+Q\^-t 



{qz + Qz) dt 



(B14) 



The integral in brackets is evaluated by introducing the 
generating function 



■E 



VBZ ^ .0 



(B15) 



One obtains (|B14p from l|B15|l by taking the second 
derivative with respect to k, and then setting k = 0. For 
non-zero i?, the integral over Q is taken easily enough. 
Completing squares, one is left with 

/ (-) e"'^«e«['^+H dt. (B16) 

VBZ Jo 

Rescaling the integration variable and taking the requi- 
site second derivative with respect to k yields 



a 2/3 



(B17) 

There is a singular contribution to the total sum from the 
term in which i? = 0. However, that term is independent 
of (f, and is, therefore, cancelled when the total expression 
with q = is subtracted. The nonvanishing contribution 
is the sum over non-zero -R's of the expression in ljB17|l . 
To this is added (|B13|I . and the result is /((j) defined 
through (|B4|) . 

Note that in all the above, the parameter /3 has not 
been fixed. It is, in fact, left undetermined, and may be 
set to speed convergence of the sums. Alternatively, it 
may be left as an internal check on the procedure, since 
the final result must be independent of (3. 




0.01 



FIG. 7; The quantity A(r) = 2A(r)l{Z'^e^ IksTU), where 
A(r) is the coefficient of fc^ in the quadratic energy of distor- 
tion of the modified fee lattice, quadratic in the displacement 
field, as a function of the aspect ratio r — Ih/lv In this plot, 
the 2-component of the wave vector q of the distortion has 
been set equal to zero. 



To see what the sum developed above yields, we split 
the wave vecto q as follows: q — zQ -\- k, where the 
projection, fc, lies in the x-y plane. Numerical results 
are consistent with an energy quadratic in displacements 
that is proportional to when Q — and k is small. 
Figure [7| shows the coefficient of k"^ in that line as a 
function of the aspect ratio r = Ih/lv Note that the 
coefficient goes through zero at a value of r that is close 
to 1. This corresponds to the "spinodal instability" ly- 
ing below the first-order transition at r « 1.1. Near the 
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spinodal the harmonic spectrum takes the form shown in 
Eq. 12. 91 The coefficient A{r) goes through a maximum as 
r increases above this value, and then decreases, tending 
to zero asymptotically as the aspect ratio becomes large. 
This latter tendency reflects the weakening of the interac- 
tion that stabilizes the counterion lattice as the distance 
between the rods grows in comparison to the distance 
between neighboring charges on a rod. In fact, it is not 
hard to demonstrate that A(r) will decay exponentially 
as r, the aspect ratio, gets large. 

To deal with terms that are third and fourth order 
in the displacement field, the procedure is the same as 
the one discussed above. One is left with lattice sums to 
perform, and the Ewald method leads to rapidly converg- 
ing numerical algorithms. The results H2.10|l - H2.13|) are 
obtained. We find for the coefficients Bi-B^ in H2.11|l . 
Bk = {Z^e^/2kBTU)Bk, where 



Bi = 0.143 

^2 = 3 X (-0.163) 

S3 = 0.010 



(B18) 
(B19) 
(B20) 



APPENDIX C: COUNTERION MELTING IS 
3DXY 

As the melting transition is approached from the 
"counterion liquid" side, we assume that there is an insta- 
bility leading to a modulated counterion charge density 
on each rod; then we allow for phase fluctuations. For 
the counterion density on each rod, we write 



p{z) = Acoii{Qz + (j){z)) 



(CI) 



Phase fluctuations disorder the charge density wave 
above and render the mean charge density on a rod sta- 
tistically uniform in the "liquid" phase. The interaction 
between counterions on a single rod will be of the form 

Making use of (|C1|) . one obtains terms of the form 

j dz j dz'V{z - z') exp [iQ{z -z') + i {<j){z) - 0(z'))] 

(C3) 

We now go to "center of mass" and "relative" co- 
ordinates. Let 



Z = 
Z = 



z + z' 
2 

z — z' 



(C4) 
(C5) 



The double integral in ljC3|l is, then, proportional to 

dZ / dZV{Z) exp {iQZ + {Z) Z + ■ ■ ■] 
= J dZ J dZV{Z)cxp[iQZ] 

X (^1 + 2i<j)' {Z)Z + ^ {2i(j>' (Z) Zf + ■ 

dZ Uq) - 2i^'{Z)v'{Q) + 2 {dp'{Z)Y v"{Q) 

(C6) 

Now, let's assume that 



v{q) (X I dz eyiYi[iqz\V {z)dz 



(C7) 



has a minimum at g = Q. Then, v'{Q) = 0, and v" {Q) > 
0. There is thus a contribution to the total energy going 
as / {d(f>{z)/dz) dz. 

Now consider the interaction between rods. One ex- 
pects that there will be terms of the form 



dzi J dz2W{zi — Z2) 

X exp [iQ (zi - Z2)+i (0i(zi) - (/i2(z2)] (C8) 



The (/>'s are subscripted to make it clear that they refer 
to different rods. If one assumes that the interaction W 
is sufficiently short-ranged, which seems to be the case 
even for unscreened Coulomb interactions because of the 
sinusoidal nature of the assumed state, one reveals the 
essence by replacing 02 (^2) by (p2{zi), and the integration 
in l|('8(l becomes 



dzi exp [i{4>i{zi) ~ 4)2izi))] 
X J dz2W{zi - Z2) exp [iQ{zi - Z2)] 
= w(Q) J dziexp [i((/)i(zi) - 02(21))] 

w(Q) cos [0i(zi)- 02(^1)] d^i (C9) 



The result of this heuristic derivation is that there are 
terms in the energy going as {d(j){z) / dz)"^ and that there 
are also terms going as cos(0i(z) — (f>2{z)), where the 
subscripts refer to near-neighbor rods. The universality 
class for the transition in this model is that of a 3DXY 
model. That the model is spatially anisotropic does not 
influence the universality class. 



